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Abstract. Certain localised displacement fluctuations in the planar honeycomb 
lattice may be identified as precursors to topological defects. We show that these 
fluctuations are among the most pronounced non-affine distortions of an elemental 
coarse graining volume of the honeycomb structure at non zero temperatures. 
We obtain the statistics of these precursor modes in the canonical ensemble, 
evaluating exactly their single point and two-point spatio-temporal distributions, 
for a lattice with harmonic nearest neighbour and next near neighbour bonds. As 
the solid is destabilised by tuning interactions, the precursor fluctuations diverge 
and correlations become long-lived and long-ranged. 


1. Introduction 

One of the most striking differences between crystalline and amorphous solids is the 
difficulty of defining analogs of lattice defects such as vacancies, interstitials and 
dislocations in the latter mm- The lack of a unique reference configuration in an 
amorphous solid precludes the use of such geometrical devices such as ideal lattice sites 
or Burgers circuits. Recent studies of the mechanical deformation of glasses do offer a 
clue, however [3]. It is now accepted that one may identify regions in the amorphous 
solid which undergo non-affine distortions mm®, i.e. displacements that cannot 
be represented as a uniform strain or rotation acting on a reference volume, £2, in the 
presence of external stresses. The dynamics of such regions has been used to describe 
irreversible plastic deformation of amorphous solids, in ways roughly analogous to 
that of dislocations in crystals m®. Continuing along this line of thinking, one may 
ask whether configurations involving topological defects in crystalline solids may be 
singled out using only measures of the local non-affine displacements [8] of atoms away 
from their ideal reference positions, without explicit reference to a Burgers circuit? In 
Ref. an affirmative answer to this question was given for the case of dislocation 
dipoles in the ideal two-dimensional (2D) triangular lattice. It was shown (1) that 
specific localised vibrational modes contribute to non-affine displacements in crystals, 
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(2) the largest contribution to thermally excited non-affine fluctuations arises from 
particle displacements that tend to introduce a dislocation dipole, i.e. localised defect 
precursor modes and (3) this contribution is separated by a large gap from the next 
largest contribution, which represents other, more complicated deformations of D. 
Here, we continue the program initiated in [9] for a lattice structure that is the 
geometric dual of the triangular lattice, namely, the planar honeycomb (PHC) crystal. 

While emphasising that our work is not necessarily specific to any particular 
material, we must mention that there are, indeed, compelling reasons to study the 
PHC lattice. This decade’s two most ground-breaking materials, with tremendous 
technological applications are Graphene nni and hexagonal Boron Nitride (h-BN) 
films El. They are known to be flat mono-layers of atoms arranged in the form of a 
two-dimensional (2D) honeycomb lattice. Apart from these, the PHC lattice appears 
in many other atomic and macro-molecular systems such as confined water [12 , 
assemblies of patchy colloids m and nano-particles M- Mechanical properties of 
lattice defects in these materials play a major role in deciding their applicability in 
various facets of the industry. Defects not only determine strength but also influence 
the topological characteristics of 2D sheets by coupling to the local curvature, leading 
to rolling up of 2D sheets into tubes and spheres or the creation of ripples on the 
surface m- Obtaining greater insight into the physics of such defects in the PHC 
lattice is therefore crucial to our understanding of these features exhibited by Graphene 
and similar materials. 

As an example, consider the Stone-Wales [T5] (SW) topological defect observed in 
sp 2 -bonded carbon materials as well as in boron nitride nanotubes and nanosheets [T71 
ITS] . A SW defect is formed due to an in-plane rotation of a carbon-carbon 
bond by 90°with respect to the bond center, resulting in bond reconstruction and 
eventual transformation of four hexagons into two heptagons and two pentagons (see 
Fig. [l]). The SW defects play a significant role in the formation of various carbon 
nanostructures, through local amorphization of Graphene sheets m- They are also 
known to influence the electronic [20j as well as transport properties |21: of Graphene, 
which can be exploited in the development of carbon-based electronic devices. 





Figure 1 . (Color online) The formation of a Stone Wales defect in graphene-like 
systems through rotation of the central basis (orange) by 90°. 

In this paper we focus on lattice deformations of the ideal PHC crystal with the 
specific objective of identifying non-affine localised modes which act as precursors 
to the formation of topological defects. A description of lattice defects in terms 
of local non-affine displacements may be useful for many reasons. Firstly, there is 
some evidence that non-affine displacements are associated with defect nucleation [22] 
and may act as a predictor for its appearance. Secondly, this is a natural choice 
for describing processes such as amorphization where a solid breaks down by the 
proliferation of lattice defects m • Lastly, this may prove to be a common language 
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for describing mechanical deformation in both crystalline and amorphous matter. 

We begin by generalising the coarse-graining procedure introduced in |8] to derive 
exact statistics of the non-affinity parameter x quantifying non-affine displacement 
fluctuations in the PHC lattice. Examining the nature of the localised vibrational 
modes responsible for producing y, we discover that one of the largest contributions 
to x arises from a mode that appears as a precursor to the SW topological defect. As 
expected, non affine fluctuations in general and the SW precursor mode in particular 
increase in amplitude as the solid is made softer. We also show that in the same 
limit, correlations between SW precursor fluctuations become stronger: the correlation 
length increases and so does the lifetime of these fluctuations. Spatial correlations 
between SW precursor modes are also strongly anisotropic and therefore encourage 
clustering, producing cascades of SW distortions of the kind described in [23]. Our 
results may be explicitly confirmed by direct use of video microscopy [24 on patchy 
colloids. These form PHC lattices, and the interaction strength may be varied using 
surface modifications m- Some of our results may also have implications for the 
mechanical behaviour and amorphization of Graphene as discussed later. 

The rest of the paper is organised as follows. In the next section we describe 
in detail the process by which we identify non-affine fluctuations of a coarse graining 
region in the PHC lattice. In order to obtain explicit and analytic results and to make 
our calculations generic, we use a simple parametrisation of the particle interactions in 
terms of harmonic nearest neighbour (NN) and next near neighbour (NNN) harmonic 
bonds. This is motivated by the fact that for small fluctuations any realistic interaction 
may be mapped onto a harmonic solid. In section 3, we describe our results for the 
single point and two point statistics of defect precursors. We end the paper with a 
summary of our main results and conclusions. 

2. Coarse graining displacement fluctuations 

At any non-zero temperature T, atomic displacements fluctuate from their reference, 
ideal lattice positions with amplitudes determined by their interactions. Given any 
sub-volume Ll within the solid, these fluctuations may be either affine or non-affine. 
The former case represents fluctuations that simply deform the entire sub-volume 
uniformly. In two dimensions, there are only four such affine transformations possible: 
a uniform change in volume, a uniaxial strain, a shear or a uniform rotation. If 
the number of particles within the sub-volume is larger than three, then there is a 
possibility of deformation modes which are non-affine, i.e. cannot be described by a 
linear combination of the four simple affine transformations of Q. The statistics of 
these non-affine displacements have previously been obtained for the one dimensional 
chain and the two dimensional triangular lattice [8]. It has also been shown that 
dominant non-affine fluctuations in the triangular lattice correspond to specific atomic 
displacements that tend to produce lattice defects f9]. We now use the ideas introduced 
in |S] and [9] to obtain the spatio-temporal statistics of non-affine fluctuations in the 
planar honeycomb (PHC) structure. 

We start by defining a simple set of interactions that can stabilise a PHC lattice. 
These cannot be standard central force potentials as these usually lead to closed packed 
solids such as the triangular lattice in two dimensions |2]. Open lattices need strongly 
orientation dependent potentials. In atomic and molecular systems this is due to 
covalent or hydrogen bonds. For example, modelling of properties of carbon based 
materials is generally performed using bond-order, short-ranged empirical potentials 
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Figure 2. (a) Coarse graining volume Q of the PHC lattice, marked in green. 
The central basis atoms (orange, labelled 1 and 2) possess 12 neighbors (blue) 
within this domain. The nearest neighbor bonds are indicated by bold lines while 
the next nearest neighbor bonds for one of the central basis particles are shown by 
dotted lines, (b) Phonon dispersion curves a;(q) for k = 1.0 and 0.4, respectively. 
The inset defines the high symmetry directions T, M and K in the Brillouin zone. 
Note that all the modes soften as k decreases. As k —> 0, the shear modulus 
vanishes and the transverse acoustic branch becomes unstable (cj(q) < 0) at 
small q. 


such as the Brenner [26] and Tersoff [27] potentials. Similar potentials are also 
used to model water [12]. In colloidal solids, one requires patchy interactions, i.e. 
heterogeneities on the surface of the colloid particles that produce angle dependent 
forces ESUHHU. Here, in order to obtain results that are sufficiently generic, we use 
a PHC network of particles connected by harmonic bonds, noting that all interactions 
are harmonic in the limit of small fluctuations. Accordingly, we assume the following 
general Hamiltonian, 


N 2 -I N 

Hharm = E + 2 ^ ^ “ U P' (1) 

i i^j 

where and u$ are the momentum and displacement (from the ideal reference) of 
particle i, m is the particle mass and is the stiffness of the harmonic bond between 
particles i and j. It is well known that a PHC lattice with only NN interactions is 
unstable [28] . To make the lattice mechanically stable, we need to add either NNN 
bonds between the particles or include an energy for bond bending. In this work, 
we choose the former path although we have checked that the two methods give 
qualitatively the same results. We therefore set Kij m K\ ~ 1 if i, j form a NN pair 
and = K 2 = n if they are next nearest neighbours. All other K; t j are set to 
zero. If we call l the lattice parameter, i.e. the length of a NN bond in the reference 
configuration, then /, \K\l 2 and \Jm/K\ set the units for length, energy and time. 
We therefore set them to unity in the following. 

For the coarse-graining analysis it will be useful to have the interaction energy 
written in a more generic form that makes it easy to switch to a representation in 
terms of Fourier modes. Looking at Fig. [2|a), which shows a portion of the PHC 
lattice, the structure may be constructed by inserting a basis of two atoms separated 
by a distance l on a triangular primitive lattice. From now on we therefore label each 
particle by the pair (i, <a), where i refers to the specific basis or unit cell, and a = 1, 2 
identifies the particle within the two-atom basis for the lattice. We call the particle 
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displacements u i a and their cartesian components u™ a with n = x,y. 

In terms of these variables, the configurational part of the Hamiltonian takes the 
generic form 

H = \ E (2) 

Here is the overall, system-wide dynamical matrix. Now define the Fourier 

transform of the particle displacements, u a (q) such that the real-space displacements 
are u i a = Vg Z f dqu a (q)e 2q Ria . Here and below the q-integral runs over the 
Brillouin zone, and the normalising factor vbz = 87r 2 /3v / 3 is its volume. The 
dynamical matrix in q-space, for the PHC lattice with NN and NNN 

interactions then has the following structure: 


Wx 

w 2 

w 3 

w 4 \ 

w 2 * 

w 5 

w 4 

W 6 

ws 

wi 

Wx 

w 2 

wi 

w 6 * 

W 2 * 

W 5 ) 


with 

3 3 v^3 

Wi = - + 3 k - 3Kcos(-q x )cos(—q y ) 

3 y/3 

W 2 = V3K,sm(-q x )sm( — q y ) 

1 1 y/3 

W 3 = - exp (iq x ) - -exp{--i q x )cos( — q y ) 

jjr / 1 • N, • / VS . 

w 4 = —iexp(--tq x )sm(—q y ) 

3 3 V^3 

w 5 = - + 3 k - 2kcos(a/3 q v ) - ncos(-q x )cos(—q y ) 

W 6 = - ^exp(-^iq x )cos(^q y ) 


(3) 


The row and column entries are arranged in the order (a,n) = (l,x), (1 ,y), (2,x), 
(2 ,y) here. The eigenvalues cj s (q) of the above dynamical matrix give the phonon 
dispersion curves for each of the s = 1... 4 phonon branches. The dispersion curves 
for k = 1.0 and k = 0.4 are plotted in Fig. [2^b) . As expected there are two acoustic 
branches where cj s (q) —>> 0 as q —>> 0 and two optical branches. As n 0 the lattice 
becomes unstable [28, to shear fluctuations and cj(q) for transverse acoustic phonons 
is imaginary at small q. 

From the dynamical matrix one can now determine the statistics of all 
displacement fluctuations. The variances of the displacements in Fourier space are 
given by 

/ m ( \ m ( !\\ a sa( C ^ a s'y( € l f ) , r\ (A \ 

{u a (q)ix 7 (q )> = 2^-- v BZ S{q + q ). (4) 


where the angular brackets indicate a thermal average and a s (q) is the eigenvector of 
q) corresponding to the s-th phonon branch. Here and in similar formulas below we 
omit an overall factor of ksT that simply scales all displacement variances. Utilizing 
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the above Fourier space variances one can then express the covariances of the real-space 
displacements as 


lipJ 1 7 / 
\ a ia +77 


>=+/ — 


dq a^(q)< 7 (q)* 
vbz ^ f ( q ) 


piq^R-io; R j7 ) 


(5) 


We now move on to the definition of coarse-grained elastic strain and non-afhne 
displacements. The green shading in Fig. [2|a) identifies the coarse-graining region Q 
we will use. It consists of a central basis (shown in orange) surrounded by a nearest 
neighbour shell of 6 x 2 = 12 atoms (shown in blue). Note that D is the smallest non¬ 
trivial neighbourhood consistent with PHC symmetry that includes all the interactions 
necessary for obtaining a stable solid. 

To define non-afhne displacements, we consider relative displacements with 
respect to either of the two atoms in the central basis. Labelling the central basis 
as i = 0 and the others i = 1,..., 6, we denote these relative displacements by A^ 7 
= u i a — uo 7 . The non-affinity then measures [5] how well these relative displacements 
are approximated by the best-fit affine deformation: 

x = min I ^2 [ A ia 7 “ D ( R ia - Ro 7 )] 2 ] (6) 

\ (ta 7 ) / 


where the minimisation is over all four elements of the two-dimensional deformation 
D. We include in the sum over triples (icry) all relative displacements across NN and 
NNN bonds involving either one of the particles in the central pair. The choice of NN 
and NNN bonds is made as these bonds carry the physical interactions. Each atom in 
the central pair has Nq = 3 + 6 = 9 such bonds, so the total number of terms in the 
sum in (J6| is 2 Nq — 1 = IT, the —1 arising from the double counting of the central 
NN bond (ia^y) = (0,1, 2). 

Expression can be concisely written [8] as 

X = min(A — Re) 2 (7) 

e 


where e is the column vector (Du, D 12 , ^ 21 , D 22 ) T formed out of the elements of the 
deformation tensor, A is a column vector collecting the relevant A^ a7 , and R is a 
[2(2 Nq - 1)] x (2 x 2) = 34 x 4 matrix with entries = S mn (R^ - Kg^). On 

carrying out the minimization in 0 one obtains 

e = QA (8) 

with Q = (R T R) -1 R T . The corresponding minimum value of the non-affinity parameter 
is given by 

x = (A - RQA) 2 = A t PA (9) 

Here P = I — RQ turns out to be a projection matrix that projects onto those 
directions of A that cannot be modelled by an affine strain. One easily checks that P 
is idempotent and hence has eigenvalues of 0 and 1, as required for a projection. 

We now proceed to obtain the statistics of x- For the particle displacements in 
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Figure 3. Log-log plot of the non-trivial eigenvalues of PCP. Note that the 
eigenvalues spectrum bunches together into band-like structures. The curves 
corresponding to the largest 3 eigenvalues are labelled by numbers. On the right, 
we show the eigen-displacements (red points and lines) corresponding to these 
eigenvalues. The reference lattice positions are shown in blue. A Stone-Wales 
defect will result from the precursor mode 2 if anharmonic interactions (neglected 
in this work) makes the NNN bonds a-b and c-d strong while breaking the NN 
bonds a-d and b-c as a consequence of this deformation. 


f2, the covariances (A™ aX A £ 7<5 ) = are given by 




£ 


/ 


dq 

VBZ l 


<C(q)a" (q)V 


iq-(Rio;—Rj 7 ) 


-a™(qK 5 (q)*e^( R -- R -) 

-a™ A (qK 7 (q)V q -( Ro *- R -) 

+ <A(q)<*(q)v q - (R -- R “ ) ' 


( 10 ) 


This Hermitian matrix with these elements, denoted C below, determines the effective 
interactions between the coarse grained variables A such that ensemble averages of 
any observable (A) defined within Q are given by 


' ' /e-i ATC - 1A dA 

In the next section, we shall use this formula to obtain coarse grained statistics of the 
non-affine displacements and of the non-affinity parameter y. 

We note two technical details: in the matrix C we use the same convention as 
already in P, Q, R and A, i.e. we include only the triples (icry) that are used in the 
definition ©• This makes the matrix C of size 34 x 34. On the other hand, with 
2 + 12 = 14 particles in Q there are only 2 x (14 — 1) = 26 independent relative 
displacement components, so C must have 8 zero eigenvalues. In expressions involving 
C -1 , the inverse should then always be understood as regularised appropriately, e.g. 
by the addition of a small constant along the diagonal of C. 
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Figure 4. P(x) and P(p) for k = 1.0 (blue) and 0.4 (red). Predictions from 
exact calculation (line) are compared with results obtained from MD simulations 
(points) of a IV = 100 X 100 PHC lattice. The system was allowed to equilibrate 
for 1 x 10 5 MD steps with a time step of 10 -3 , after which configurations were 
collected over a maximum of 8 X 10 5 MD steps at intervals of 100 time steps. 


3. Results: Statistics of non-affine displacements 

The non-trivial eigenvalues Gj > 0 of PCP correspond to non-affine displacements. 
Accordingly we have from for the mean non-affinity 

(x) = Tr(PCP) = J2 a r 0 1 ) 

3 


and for its variance, 


(X 2 ) - (X) 2 =2Tr(PCP) 2 = 2]T4 (12) 

3 

To obtain the statistics of y, we need to evaluate the Brillouin zone integrals that define 
C numerically. We use a 64 point Gauss-Legendre quadrature routine to perform all 
such numerical integrations over the Brillouin zone in this paper. The computed PCP 
is a 34 x 34 matrix with 34 eigenvalues. As explained above this has 8 eigenvalues that 
are automatically zero, and a further 4 that represent the four affine deformations in 
2D, namely volumetric, axial, shear and rotation; recall that P projects into the space 
orthogonal to these deformations. The remaining 22 non-trivial eigenvalues relate to 
distinct non-afhne modes and are plotted in Fig.[3]as a function of k. All our numerical 
results are given for a dimensionless temperature of k^,T = 0.01, where a harmonic 
approximation should be reliable provided that anharmonic interaction coefficients are 
of order unity in our dimensionless units. 

As k 0, most of the eigenvalues diverge as ^ . This behavior is expected 

since the q-integral defining C involves the inverse of co^(q), which vanishes for the 
soft phonon branch as n —> 0. Remarkably, the eigenvalue spectrum contains large 
gaps with the three largest eigenvalues being separated from the others by almost an 
order of magnitude unless k is close to unity. The eigenvectors corresponding to these 
three eigenvalues are also shown in Fig. [3J Out of these, two deformations correspond 
to breathing modes, where the central basis is shifted relative to the others in Q in the 
horizontal and vertical directions. The eigenmode corresponding to a 2 is interesting: it 
represents a rotation of the central basis relative to the reference positions. A rotation 
of ±90° together with a reconnection of the bonds would nucleate a SW defect. This 
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Figure 5. C X (R, 0) and C p ( R, 0) for k = 1.0 and 0.4. Note that while C X (R, 0) 
is isotropic, Cp(R, 0) shows prominent quadrupolar anisotropy. Both correlations 
increase with decreasing k. 


is, then, one of the non-affine precursor fluctuations that must precede nucleation 
of SW defects. Our calculations show that such precursor fluctuations are the most 
prominent (and softest ) non-affine distortions that tend to change the lattice topology. 


3.1. Probability distribution functions 

To obtain the probability distribution of y, P(y), we begin with the characteristic 
function [8], 


m = J d X P(x)e ikx = (e ikx ) 


J e -lA T (C- 1 -2i/cP)A dA 

/ e _ ^ ATC_1A dA 


(13) 


Here, k is the variable conjugate to y. By carrying out the Gaussian integrals, we 
obtain 

0(Jfe) = |I-2i/cPCP|- 1/2 . (14) 

The PCP matrix above is Hermitian and can be diagonalized by a unitary 
transformation. Hence, writing the determinant in terms of products of its eigenvalues 
we can also write 

(j){k) = W&j 2 where f>j = 1 — 2 ik<jj (15) 

3 


A numerical inverse Fourier transform of <fi(k) then gives P(y). At a more intuitive 
level, one observes from the form of ( fl5| ) that P(x) i s the distribution of the sum of 
the squares of uncorrelated Gaussian random variables with variances (jj. 

To determine the probability distribution of the SW mode we proceed similarly. 
Let p be the projection of any A in 11 along the SW eigenmode, i.e. p = bJ^A 
where bsw is the eigenvector of PCP corresponding to the SW eigenvalue. As A has 
a zero mean Gaussian distribution with covariance matrix C, p is also Gaussian with 
distribution 


p( P ) = 


i 


2tt bLCb 


: exp 


u sw 


£ V 

2b^|y C b sw J 


(16) 


The probability distributions for x an d P are shown in Fig. [4ja) and (b), where 
we have compared the analytical results, obtained by Fourier transforming the 
corresponding characteristic functions, with data obtained from molecular dynamics 
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Figure 6. (a) C x (0, t) and (b) C p (0,t) for k = 1.0 (blue) and 0.4 (red). Points 
are from simulations of a 500 X 500 PHC lattice while the lines are from our exact 
results. The simulation data was obtained from 2 x 10 3 MD configurations and 
further averaged over shifts of the origin of the time axis. 


(MD) simulations using the parallelized MD package LAMMPS [29]. Note that both 
(x) and ( p 2 ) increases as k. It is straightforward to show that Pb^vp = bsw, he. bsw 
belongs to the non-trivial eigenspace of P, implying immediately that (p 2 ) = ~ ft -1 * 

The leading behaviour of (y) as k 0 is the same. 


3.2. Correlation functions 

In order to calculate the spatio-temporal correlations of the non-affinity x and the 
SW projection p, we need to consider simultaneously displacement differences in two 
neighborhoods fi and Q centered at Ro and Ro, at times t and t' respectively. 
The non-affinity and SW projection are defined as y(Ro,£) = A T (t)PA(t) and 
p(Ro,t) = b^^A(t) for the coarse graining volume Q, with particle displacements 
within this region at a time t being represented by the vector A (t). For the region 
9 at a time t\ we correspondingly have the quantities x(Ro,t') = A T (t')PA(t') and 
p(Ro,t 7 ) = b^ w A (£'). The displacement correlations in this case are defined as 

<A(t)A T (t')) = C (17) 


where the matrix C has elements given by the expression 

<JqK4q)v q - (R -- ft ^ 


s-rmn 


-C(q»(q)*e iq ' (Ri “- ft » !) 

- <«S(qW 1 (q)*e iq ' (Roi " fi ” ) 

+ aS(qK { (q)*e iq,(R “ i - a " <) 


cos[a; s (q)(t / - t)} 

w I(q) 


(18) 


The correlation between %(R 0 ,t) and x(R-oT) can be calculated using Wick’s 
theorem as [8], 


C x (R 0 ,t,RoV) = (x(Ro,4x(Rj, 0) - (x) 2 

= 2 Tr(PCP)(PCP) T 

= 2 

3 
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Figure 7. C P (R, t) calculated at different times t = 0,5,10 and 15, for the 
harmonic honeycomb lattice with k = 0.4. Note the development of strong spatial 
anisotropy with time. 


where the <j| denote the eigenvalues of the matrix (PCP)(PCP) T . Similarly, the 
correlation of p( Ro,£) and p(Ro,t') is given by, 

C p { Ro,t,R 0 ,t') = (b T sw A(t)A T (t')b S w) 

= t>5vv(^(^)^ T (^ / ))^ ) *s'w (20) 

= bsw 

In a homogeneous solid in equilibrium, these correlation functions are functions only 
of the relative coordinates Ro — Ro and times t — t'. From now on we will simply 
denote these by R and t respectively. 

We first focus on the spatial correlations at equal time. These are depicted in 
Fig. [5} where we plot C X (R, 0) and C P (R, 0) for two values of k. The y correlation 
are short ranged, isotropic functions qualitatively similar to those for the triangular 
lattice UM and increase in magnitude as n decreases. The SW precursor correlations 
are also short ranged but show prominent quadrupolar anisotropy. The range and the 
anisotropy increases grows as the PHC lattice softens with decreasing n. Non-afhne 
distortions are localised in space being composed of many incommensurate lattice 
distortions. The SW precursor mode, can condense into SW defects when an-harmonic 
terms in the interactions become prominent in the small n regime. 

In Fig. [6] we turn to the local correlations as a function of time. The figure 
shows the normalized correlations, C x ( 0,£) and C p ( 0,£). Our results are compared 
with those obtained from molecular dynamics simulations of a 500 x 500 site harmonic 
honeycomb net. As is evident from the figure, the correlations decay with time showing 
that they are transient. The decay occurs in an oscillatory fashion due to the presence 
of a large number of mutually incommensurate frequencies. The mean frequency of 
these oscillations decreases as k 0, indicating larger relaxation times of the induced 
deformation. 

Finally, to demonstrate the behavior of the correlations over both space and time, 
we plot C P (R, t) in Figj7| for k = 0.4, at four specific values of t. Interestingly, the 
spatial anisotropy of the correlations increases with time even as their amplitude 
decays. The behaviour of C x ( R, £) is similar to that observed in the triangular 
lattice 

4. Discussion and Conclusions 

In this paper we have carried out a detailed study of the displacement fluctuations 
of a coarse-graining volume, Q, within the PHC crystal. At non-zero temperatures, 
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thermal fluctuations produce random affine local deformations (such as scaling, shear 
or rotation), but also distortions of Ll that are impossible to represent as an affine 
transformation. We have obtained the full single point and two-point statistics of 
these non-affine fluctuations in the harmonic limit. The approach to the limit where 
the lattice is unstable causes a divergence of the mean non-affinity, and of the spatio- 
temporal scales of the correlation functions. We have shown that the vibrational 
modes that contribute to non-afhne fluctuations are given by the eigenvectors of the 
matrix PCP. This is constructed from displacement correlations C and the operator P 
that projects all atomic displacements onto the non-affine subspace. The eigenvalues of 
PCP are seen to group together into three well defined clusters with large gaps between 
them. As the lattice softens, most eigenvalues diverge but the clusters remain distinct. 
We have shown that one of the three largest eigenvalues involves a rotation of the basis 
unit relative to the rest of the atoms in Ll and therefore would tend to nucleate a SW 
defect in a system with a more realistic interaction potential. We have obtained the 
statistics and correlations of these localised SW precursor modes in the PHC lattice. 

Most work on lattice defects in the PHC crystal have centered on Graphene for 
obvious reasons and have involved mainly the calculation of defect formation energies 
and energy barriers at zero temperature. Within our model, we estimate that a 
value of n = 0.4 — 0.2 reproduces the phonon spectrum of Graphene reasonably well. 
Comparing the velocities for acoustic phonons for k = 0.4 with the actual values in 
Graphene, we also estimate that our dimensionless temperature corresponds to a real 
temperature of about 290K. Our predictions for the spatial and temporal correlations 
for x and the SW precursors should be accurate for these parameters. SW defects 
in Graphene typically have a very high activation energy 5.39 eV) as determined 
by ab-initio calculations m ■ This can be lowered by extrinsic mechanisms such as 
hydrogen adsorption, but only to about 2.54 eV m- This would make SW defects 
unlikely in pristine Graphene at temperatures which are relevant. Our calculations 
suggest that this picture may need to be modified. 

Entropic contributions to defect stability, arising especially from lattice vibrations 
that have, to the best of our knowledge, received relatively less attention, appear to be 
at least equally important. This is not particularly surprising, in retrospect, because 
in 2D the loosely packed and open PHC lattice is expected to be strongly influenced 
by thermal fluctuations ^ 32]. 

Our calculations also show that SW defects are highly correlated in space and 
time due to elastic interactions. Interpreting the space-time correlation functions 
C P (R, t) as a linear response function 0 0 |33j, one can argue that an isolated SW 
defect at any point R would induce other SW defects at neighbouring lattice points. 
The spatial distribution of these induced SW defects would be anisotropic and follow 
C p ( R, £) as shown in Fig. [5j Similarly, Fig. [7| shows that the anisotropy in the SW 
precursor correlations increases with time. This may be interpreted as implying that 
SW defects are nucleated in cascades [23], which enables these defects to proliferate 
during amorphization m — at which point the harmonic theory breaks down. 

How does external stress X) influence x and the SW precursors? To answer this 
question one only needs to include [8j|9] the term T, t J2iLi e(R$) in the Hamiltonian 
0. where e(R$) is the best-fit local strain, QA, at particle i. To lowest order in X, the 
variation of x is then (x)s = (x)s=o + ^ 5] T QC T PCQ T 5], where the dependence on 
i enters in C via the relative displacement between the location where stress is applied, 
R^, and the location where we are measuring x, say Ro- The probability P(p) remains 
Gaussian with the same variance but now the db symmetry of p is broken and P(p) 
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is shifted with a mean (p) = ^bc^CG^X) 7 ^ 0 . Stress, X, also does not affect the 
space-time correlation functions of p. The density of SW defects will be modified, 
though, by curvature terms m arising from out of plane fluctuations that have not 
been included in our calculation. Curvature then acts as an external field coupling 
to x and p. The effect of external fields on x has been studied in the context of the 
triangular lattice [9] and should lead to similar results here. 

Since microscopic and instantaneous atomic coordinates are indeed difficult to 
obtain for systems such as Graphene and h-BN, we believe that our predictions may 
be tested in future using video microscopic data m on patchy colloids mu Eg that 
form a PHC lattice. Similar experiments have been successfully carried out in simpler 
colloids [24j [34l E2] so it would be interesting to see results for the more complex PHC 
structure. 

In future we would like to investigate the three dimensional FCC lattice where a 
similar association of non-affine displacements with defects has been observed [22 . In 
three dimensional close packed solids, the elementary topological defect is a small 
dislocation loop which may be difficult to identify using a Burgers circuit for a 
fluctuating crystal at T 7 ^ 0. However, local non-affinity should remain an easily 
recognised marker. Analysis of non-affine modes in progressively disordered solids 
should also elucidate the connection, if any, between topological defects in crystals 
and non-affine droplet fluctuations or shear transformation zones in amorphous 

matter urn- 
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